CETN-I-43 

6/87 


(?o44fal  Sttyittewiitty 

t7ec6*ucal  'Hote 

SPECTRAL  ANALYSIS  OF  DIGITAL  WAVE  DATA 
COMPUTER  PROGRAM:  SPECTRUM 


Purpose :  This  CETN  introduces  an  interactive  program  for  spectral  analysis  of 
water  elevation  data  capable  of  fast  execution  on  microcomputers. 


Capabilities :  Computer  program  SPECTRUM  takes  a  time  series  of  evenly  spaced 

water  elevation  data  and  breaks  it  down  into  its  surface  variance  spectral 
density.  Under  linear  wave  theory,  the  surface  variance  spectral  density  can 
be  converted  to  spectral  energy  density  upon  multiplication  by  pg  -  where  p 
is  the  density  of  water  and  g  is  acceleration  due  to  gravity.  SPECTRUM  is 
written  in  FORTRAN  77  and  can  handle  a  time  series  of  up  to  4096  points. 
SPECTRUM  is  very  efficient  in  terms  of  memory  usage  and  thus  is  attractive  for 
use  on  microcomputers. 


SPECTRUM  will  test  for  a  trend  in  the  time  series  and  optionally  delete  it 
upon  consent  of  the  user.  SPECTRUM  also  produces  90  percent  confidence 
intervals  for  the  actual  surface  variance  spectral  density  values.  SPECTRUM 
calculates  the  zeroth  through  fourth  moments  of  the  spectrum  and  several 
useful  wave  and  spectral  parameters  including  significant  wave  height  and  peak 
period.  SPECTRUM  also  allows  the  user  a  choice  of  data  windows  or  tapers. 

The  interactive  input  has  optional  default  settings  for  most  of  the  input 
data. 


Methods :  The  data  analysis  routine  is  performed  using  the  Welch  method 
(Welch,  1967)  of  spectral  estimation  which  may  be  briefly  summarized  as 
follows: 

(1)  Segment  the  N  point  time  series  (Yn)  n  =  0,1 . N-l  into  K 

half -overlapping  segments  of  L  points  each.  This  implies  that 


K  =  (N  -  L/2)  /  (L/2)  and  we  define  the  jth  segment  as 

Z.  =  Y  t  /o  where 

jn  n+(j-l)  L/2 

j  =  1,2  ....  ,K 
n  =  0,1,  .  .  .  .  , L— 1 
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(2)  Subtract  out  the  mean  of  each  segment 

X  =  Z.  -  Z.  where 

jn  jn  j 

j  1,2,  •  •  *  •  ,  K 
n  —  0,1,  •  •  •  •  , L- 1 


Z .  =  L E,!  Z.  /L 
J  n=0  j  n 

(3)  Compute  G.(mAf),  the  discrete  Fourier  transform  of  the  jth  segment  at 
frequency  ^mAf 

G.(mAf)  =  At  X.  W  e_2lTlmn/L 

j  n=0  j  n  n 

j  =  1,2  .  »  •  .  ,K 


m  =  0,1,  .  .  .  .  ,L/2-l 
i2=  -1 


Af  =  1/  (NAt) 

At  =  the  recording  time  between  consecutive  data  points 
W  =  the  data  window  for  the  nth  point  of  the  segment 


(4)  Compute  the^estimated  surface  variance  spectral  density  from  the  jth 
segment  -  Sj(mAf) 

s\  (mAf )  =  2Cj  [r  Gj  (raAf  )|  2/LAt  where 

j  =  1,2,  .  .  .  .  ,  K 
m  =  0 , 1 ,  .  .  .  .  ,  L/ 2-1 


C.  = 
J 


the  complex  modulus 


(length)  of  Gj(mAf) 


(5)  Average  the  surface  variance  spectral  density  from  each  segment  to 
obtain  a  final  estimate  which  we  shall  denote  by  P(mAf). 


P ( m Af  )  =  .  E,  S.  (mAf)/K 
J  =  1  J 

The  principle  advantage  of  the  Welch  method  in  comparison  to  other  non¬ 
overlapping  segmenting  techniques  is  that  it  reduces  the  width  of  the 
confidence  interval  for  the  actual  surface  variance  spectral  density  values. 
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Additional  details  pertaining  to  the  mathematics  of  spectrum  and  time  series 
analysis  may  be  found  in  Bendat  and  Piersol  (1971),  Borgman  (1972),  Borgman 
(1973),  Lund  (1986),  and  Welch  (1967). 

The  discrete  Fourier  transform  is  performed  using  the  fast  Fourier  transform 
algorithm  which  is  an  extremely  efficient  computational  scheme  for 
G  (mAf).  The  fast  Fourier  transform  routine  in  the  program  requires  that 
thl  time  series  length  N  be  a  power  of  2  -  if  not,  the  time  series  length  is 
truncated  to  the  largest  power  of  2  less  than  N.  A  90  percent  confidence 
interval  for  the  actual  surf aceAvariance  spectral  density  value  P(mAf)  is 
calculated  using  the  value  of  P(mAf)  and  the  chi-square  probability  law.  The 
user  may  easily  adjust  SPECTRUM  for  percent  values  other  than  90  percent  by 
changing  the  upper  tail  Z  ,  value  from  the  standard  normal  distribution  at 
runtime.  The  zeroth  througn  the  fourth  moments  of  the  spectrum  are  calculated 
by  the  program.  The  rath  moment,  Uffl»  is  defined  as  follows: 

U  =  I  f^f )  df  m  =  0,1,2,  ...  . 

ra  Jo 

The  spectral  parameters  are  numbers  which  quantify  the  shape  of  the  wave 
spectrum.  They  are  usually  functions  of  the  moments.  The  spectral  parameters 
calculated  in  SPECTRUM  are  listed  and  defined  as  follows: 

The  spectral  width  parameter  =  e  =yj 1-  (U^^/U^U^) 

The  spectral  narrowness  parameter  =  v  = 

The  spectral  peakedness  parameter  =  Qp  =  (2/U„^)  f ”f  P(f j  df 

JO 

The  significant  wave  height  =  Hg  =  4V"U^~ 


These  parameters  are  discussed  in  Rye  (1979).  Note  that  Hg  by  this  definition 
is  not  necessarily  equivalent  to  the  average  height  of  the  one-third  highest 
waves,  particularly  in  shallow  water.  The  integrals  in  the  moments  and 
equations  defining  the  spectral  parameters  are  calculated  using  the 
trapezoidal  rule  as  the  numerical  integration  scheme  and  ^(raAf)  as  an 
approximation  for  P(mAf). 

The  peak  frequency  is  calculated  by  two  different  methods  and  is  defined  as 
the  frequency  at  which  the  maximum  surface  variance  spectral  density  occurs. 
The  first  method  searches  the  array  of  variance  density  estimates  for  its 
maximum  value;  the  second  method  averages  all  frequencies  which  have  variance 
densities  that  exceed  0.6  times  the  peak  frequency  calculated  by  the  first 
method.  The  peak  period  is  computed  as  the  reciprocal  of  the  peak  frequency 
computed  in  the  first  method. 
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The  data  windows  are  used  in  the  reduction  of  side  lobe  leakage.  This  helps 
keep  variance  from  "leaking"  into  frequencies  where  it  does  not  belong. 
SPECTRUM  gives  the  user  a  choice  of  three  data  windows: 

1)  The  rectangular  boxcar  window 

W  =  1  0  <  n  <  L-l 


2)  The  10  percent  cosine  bell  window 

/  1  / 2  [  1  -  cos(  10xn/L)  ]  0  _<  n  _<  L10 

W  -  <  1  L10+1  <  n  <  L-L10-2 

„ 

V  1/2  [  1  -  cos(  10ir(L-n-l  )/L)  ]  L-L10-1  <  n  <  L-l 

where  L10  =  [L/10]  =  the  greatest  integer  less  than  or  equal  to  L  divided  by 

10. 

3)  The  Hamming  window 

W  =  .54  -  .46  cos  (2xn/L)  0  <  n  <  L-l 

n  —  — 


Input:  A  title  for  each  run  is  entered  followed  by  the  number  of  data  points 
N,  the  segment  length  L,  the  data  window  choice,  and  the  sampling  frequency 
fs  (At  =  1/fs).  Default  values  for  the  segment  length  and  data  window  type 
are  provided  by  the  program.  The  default  data  window  is  the  10  percent  cosine 
bell;  while  the  default  segment  length  is  chosen  such  that  the  number  of 
segments  is  seven.  These  values  can  be  easily  changed  by  the  user  at 
runtime.  The  sampling  frequency  can  be  interpreted  as  the  number  of  data 
points  recorded  per  second. 

Next,  SPECTRUM  asks  for  lower  and  upper  cutoffs  for  the  computation  of 
spectral  moments  and  parameters  (Rye,  1979).  Cutoffs  should  be  in  terms  of 
factors  times  the  peak  frequency.  Default  cutoffs  are  0.3  and  3.0  times  the 
peak  frequency.  Finally,  the  units  of  water  elevation  data  and  the  name  of 
the  data  file  containing  the  water  elevation  data  are  input.  The  user  should 
make  sure  that  the  format  statement  in  the  program  reads  the  input  time  series 
in  the  correct  format.  This  should  be  done  by  the  user  at  runtime.  The 
default  format  statement  reads  the  data  in  a  6X,25F3.0  formatted  data  file. 
This  means  the  data  file  will  contain  six  blank  spaces  at  the  beginning  of 
each  row  of  data  and  will  be  followed  by  25  three  digit  floating  point  numbers 
with  no  blanks  between.  If  the  data  file  is  not  in  this  format,  the  user  will 
be  allowed  to  enter  the  format  appropriate  for  the  specific  data  file.  If  a 
trend  in  the  data  is  detected,  the  user  is  informed  and  has  the  option  of 
deleting  it  or  leaving  it  stand  as  is;  if  no  trend  is  found,  the  program 
continues  execution  without  interuption.  All  input  time  series  should  be 
reviewed  and  edited  for  quality  before  use  with  this  program.  Input  data  must 
be  water  level  measurements.  Erroneous  results  may  be  produced  if  this 
program  is  used  to  process  pressure  or  current  time  series. 
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Output :  The  title  of  the  run  is  printed  as  a  heading  for  the  program 

output.  This  is  followed  by  a  table  of  frequency  verses  surface  variance 
spectral  density.  This  table  also  includes  90  percent  confidence  intervals 
for  the  actual  surface  variance  spectral  density  values.  The  confidence 
intervals  may  be  interpreted  as  follows  -  actual  surface  variance  spectral 
density  falls  between  the  lower  and  upper  bounds  with  90  percent  confidence. 
The  table  is  followed  by  output  of  the  zeroth  through  fourth  spectral 
moments.  Next,  the  degrees  of  freedom  for  the  chi-square  distribution  and  the 
recording  interval  between  consecutive  data  points  (At)  are  printed.  Then 
the  spectral  parameters  are  output.  The  significant  wave  height,  spectral 
peakedness  parameter,  spectral  width  parameter,  and  spectral  narrowness 
parameter  are  all  displayed  by  the  program.  Finally,  peak  frequency  and  its 
reciprocal,  peak  period,  are  displayed  by  SPECTRUM. 


Sample  Problem:  A  twenty-minute  time  series  from  a  Waverider  buoy  in  Lake 
Michigan  (South  Haven,  Michigan)  is  analyzed  in  the  following  sample  execution 
of  SPECTRUM. 

8>SPECTRUM 


*  SPECTRUM  IS  A  PROGRAM  WHICH  BREAKS  DOWN  A  TIME 

*  SERIES  OF  WATER  ELEVATION  DATA  INTO  ITS  SURFACE 

*  VARIANCE  DENSITY  SPECTRUM.  SPECTRUM  IS  ABLE  TO 

*  CHECK  FOR  AND  OPTIONALLY  ELIMINATE  TRENDS  IN  THE 

*  OATA.  SPECTRUM  USES  THE  WELCH  METHOD  OF  HALF - 

*  OVERLAPPING  SEGMENTS  AND  WILL  PRODUCE  CONFIDENCE 

*  INTERVALS  FOR  THE  ACTUAL  SPECTRAL  VALUES.  THE 

*  ZEROTH  THROUGH  THE  FOURTH  SPECTRAL  MOMENTS  AND 

*  SEVERAL  USEFUL  SPECTRAL  PARAMETERS  ARE  INCLUDED 

*  IN  THE  PROGRAM  OUTPUT. 


INPUT  THE  SITE  OF  DATA  COLLECTION  OR  THE  TITLE  OF  THIS  RUN. 
TRIAL  RUN 

HOW  MANY  DATA  POINTS  DO  YOU  HAVE  IN  YOUR  TIME  SERIES? 

4096 

DO  YOU  WISH  TO  ENTER  A  SEGMENT  LENGTH  (Y/N)? 

IF  NOT,  ONE  WILL  BE  SELECTED  FOR  YOU  SUCH  THAT 
THE  NUMBER  OF  SEGMENTS  EQUALS  7. 

N 

DO  YOU  WISH  TO  ENTER  A  WINDOW  TYPE  (Y/N)?  IF  NOT, 

THE  10%  COSINE  BELL  WINDOW  WILL  BE  CHOSEN  FOR  YOU. 

Y 

RECTANGULAR  BOXCAR  WINDOW  «  1 

10  PERCENT  COSINE  BELL  WINDOW  *  2 
HAMM  I NG  W I NDOW  -  3 

SELECT  1 ,  2 ,  OR  3 
1 

HOW  OFTEN  ARE  THE  DATA  POINTS  SAMPLED? 

USE  UNITS  OF  RECORDINGS  PER  SECOND  (HZ). 

4 

DO  YOU  WISH  TO  ENTER  HIGH  AND  LOW  FREQUENCY  CUTOFFS 
FOR  THE  COMPUTATION  OF  THE  SPECTRAL  MOMENTS  AND 
PARAMETERS(Y/N)? 

DEFAULT  VALUES  ARE  0.3  AND  3.0  TIMES  THE  PEAK  FREQUENCY 
N 

WHAT  ARE  THE  UNITS  OF  YOUR  WATER  ELEVATION  DATA? 

FEET . 1 

I  NCHES . 2 

CENTIMETERS.  .  .  .3 
SELECT  1 ,  2 ,  OR  3 . 

3 

WHAT  IS  THE  NAME  OF  YOUR  DATA  FILE? 

BE  SURE  TO  INCLUDE  THE  EXTENSION. 


C  :  SHAV 1 720 . DAT 
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THE  DEFAULT  FORMAT  IS  (6X.25F3.0) 

IS  YOUR  DATA  IN  THIS  FORMAT  (Y  OR  N)  Y 

YOUR  DATA  INDICATES  A  LINEAR  TREND  WITH  AT  LEAST  95 
PERCENT  CONFIDENCE.  DO  YOU  WANT  TO  REMOVE  IT  (Y/N)? 

Y 

TRIAL  RUN 

SPECTRAL  DECOMPOSITION  OF  THE  TIME  SERIES: 

90  PERCENT  CONFIDENCE  INTERVALS 
FOR  THE  ACTUAL  SPECTRAL  VALUES 

SPECTRAL 


FREQ ( HZ ) 

DENS  1 TY 

CM*  2  /HZ 

LOWER 

UPPER 

.047 

583 . 369 

317.320 

1501 .989 

.051 

962.267 

523.418 

2477.529 

.055 

203 . 701 

1 10 . 802 

524 . 464 

.059 

51.610 

28.073 

132 . 879 

.063 

30.498 

16.589 

78.522 

.066 

43.365 

23.588 

1 1 1 .650 

.070 

20.708 

1 1 .264 

53.318 

.074 

81 .608 

44.390 

210. 1  15 

.078 

45.254 

24.615 

116.514 

.082 

25.424 

13 . 829 

65 . 459 

.086 

28.321 

15.405 

72.917 

.090 

36.365 

19.781 

93 . 629 

.094 

32.069 

17.443 

82 . 566 

.098 

152.352 

82.871 

392.257 

.  102 

85.339 

46.420 

219.721 

.  105 

84.036 

45.711 

2  16. 365 

.  109 

95.01 1 

51.681 

244.624 

.113 

158.233 

86.070 

407 . 398 

.117 

349.956 

190.356 

901 .024 

.121 

480.614 

261 .426 

1237.426 

.  125 

1080.516 

587.739 

*  2781 . 982 

.  1  29 

1840.780 

1001 .280 

4739 . 4 1 9 

.  133 

3501 .738 

1904.746 

9015.853 

.  137 

4975.678 

2706.486 

12810. 780 

.14  1 

1 1350. 140 

6173. 829 

29222 . 960 

.  145 

10564 . 840 

5746 . 672 

27201 .080 

.  148 

10258.390 

5579.980 

26412.070 

.  1  52 

12410. 590 

6750.656 

3 1 953 . 300 

.  1  56 

17176.920 

9343 . 267 

44225.060 

.  160 

25040.200 

13620.440 

64470.480 

.  164 

24908.960 

13549.060 

64132 . 600 

.  1  68 

10246 . 520 

5573.522 

26381 . 500 

.  172 

14238.730 

7745.057 

36660. 160 

.  176 

8955.002 

4871 .01 1 

23056. 260 

.  180 

7385 . 563 

4017.325 

19015. 460 

.  184 

5700.385 

3100.685 

14676.660 

.  188 

3712.231 

2019 . 242 

9557 . 804 

.191 

3798.051 

2065 . 923 

9778 . 763 

.  195 

2936.927 

1597.521 

7561 . 646 

.  199 

3033 .121 

1649.845 

7809 .314 

.  203 

2535 . 227 

1379.019 

6527 . 396 

.207 

3630.544 

1974.809 

9347.486 

.211 

2704.000 

1470.822 

6961 . 933 

.215 

3080.744 

1675 . 749 

7931 . 927 

.219 

975.079 

530.388 

2510.517 

.  223 

2411. 549 

1311. 745 

6208 . 966 

.  227 

1 527 . 485 

830 . 865 

3932 . 785 

.230 

1072 . 366 

583 . 306 

2760 . 999 

.  234 

1 679.959 

913 . 802 

4325 . 356 

.  238 

2118. 298 

1152. 233 

5453 . 938 

.  242 

1 992 . 581 

1083 . 851 

5130.257 

.  246 

2086 . 398 

1134. 882 

5371 . 806 

.  250 

791  .320 

430 . 433 

2037 . 396 

.254 

638 . 298 

347 . 198 

1643 .412 

.258 

2046.385 

1113.117 

5268 . 785 

.262 

1157.318 

629 .515 

2979 . 722 

.  266 

1 480.093 

805.086 

3810 . 765 

.270 

1685.499 

916.815 

4339.619 

.273 

1039.320 

565 .331 

2675 .916 

.277 

619.954 

337 . 220 

1  596 .  1 84 

.281 

1036.827 

563.975 

2669.496 
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.  285 

1399 . 93 1 

761 .482 

3604 . 372 

.  289 

1691 .631 

920. 151 

4355 . 409 

.  293 

393.281 

213 . 922 

1012.572 

.  297 

268.335 

145 . 959 

690.877 

.  301 

723.807 

393 .710 

1863 .570 

.  305 

513.601 

279.370 

1322.357 

.  309 

774 . 767 

421 .429 

1  994 . 777 

.313 

388.513 

21 1 .329 

1000.296 

.316 

541 .822 

294.720 

1  395.018 

.320 

430. 128 

233 . 965 

1  107 . 44 1 

.  324 

598. 140 

325.354 

1540.020 

.328 

301 .869 

164.200 

777.216 

.332 

588 . 470 

320.094 

1515.120 

.  336 

636.964 

346 . 472 

1639.978 

.  340 

619.946 

337.215 

1596. 162 

.  344 

245.769 

1  33 . 684 

632.775 

.  348 

412.251 

224.241 

1061 .415 

.352 

399.952 

217.551 

1029.749 

.355 

340.992 

185.480 

877 . 945 

'  .359 

219.530 

119.412 

565.218 

.363 

402.318 

218 . 838 

1035.839 

.367 

524.848 

285 . 487 

1351 .314 

.371 

410.288 

223 . 173 

1056 . 361 

.375 

200.606 

109. 1 18 

516.495 

.379 

187.323 

101 . 893 

482 . 298 

.  383 

168.245 

91.516 

433 . 177 

.387 

184.571 

100.396 

475.212 

.391 

137.873 

74 . 995 

354 . 978 

.  395 

217.429 

118. 269 

559.810 

.398 

228.819 

124.465 

589. 136 

.  402 

166.166 

90.385 

427.825 

.406 

1 10.602 

60 . 161 

284.764 

.410 

94.535 

51.421 

243.396 

.414 

102.065 

55.517 

262.784 

.418 

250 . 974 

136.516 

646. 179 

.422 

157.331 

85.579 

405.078 

.426 

170.032 

92 . 488 

437.778 

.  430 

131 .411 

71 . 480 

338.342 

.434 

93.565 

50.894 

240.900 

.438 

115.813 

62 . 996 

298 .181 

.  441 

81 .882 

44 . 539 

210.819 

.  445 

90.164 

49.044 

232. 144 

.449 

148.577 

80.818 

382 . 539 

.453 

129.928 

70 . 674 

334 . 523 

.457 

43.477 

23.649 

111. 940 

.  461 

97.676 

53.130 

251 . 484 

.  465 

1 10 . 491 

60 . 101 

284 . 479 

.469 

154.909 

84.262 

398 . 841 

.473 

47.188 

25 . 668 

121 .495 

.  477 

45.413 

24 . 702 

116. 924 

.  480 

31 .338 

17.046 

80 . 685 

SPECTRAL  MOMENTS  ARE 

AS  FOLLOW 

THE  0TH 

MOMENT  OF  THE  SPECTRUM 

IS 

934  . 

347 

THE  1 TH 

MOMENT  OF  THE  SPECTRUM 

1  S 

171  . 

784 

THE  2TH 

MOMENT  OF  THE  SPECTRUM 

IS 

34  . 

466 

THE  3TH 

MOMENT  OF  THE  SPECTRUM 

IS 

7  . 

726 

THE  4TH 

MOMENT  OF  THE  SPECTRUM 

1  s 

1  . 

968 

THE  DEGREES  OF  FREEDOM  ARE 

9 

.  800 

THE  RECORDING  INTERVAL  IS 

.250  SECONDS 

SIGNIFICANT  WAVE  HEIGHT  IS 

1  22 

.268  CM 

spectral  peakedness  parameter  i 

s 

1 

.965 

THE  SPECTRAL  WIDTH  PARAMETER  IS 

.  595 

SPECTRAL  NARROWNESS  PARAMETER  IS 

.  302 

THE  PEAK  FREQUENCY  IS 

.160  HZ 

THE  DELFT  PEAK  FREQUENCY  IS 

.160  HZ 

THE  PEAK  PERIOD  IS 

6 

.244  SECONDS 

Stop  -  Program  terminated. 
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Program  Availability:  SPECTRUM  is  available  in  Microsoft  FORTRAN  for  the  IBM- 
PC  on  a  5-1/4  in.  diskette  or  as  a  printed  program  listing  and  may  be  obtained 
from  Ms.  Gloria  J.  Naylor  at  (601)  634—2581  or  FTS:  542—2581,  Engineering 
Computer  Programs  Library  Section,  U.S.  Army  Engineer  Waterways  Experiment 
Station,  P.0.  Box  631,  Vicksburg,  MS  39180-0631.  Questions  concerning 
applications  of  SPECTRUM  can  be  directed  to  (601)  634-2012  or  FTS:  542-2012 
the  Coastal  Design  Branch,  Coastal  Engineering  Research  Center,  U.S.  Army 
Engineer  Waterways  Experiment  Station,  P.0.  Box  631,  Vicksburg,  MS  39180- 
0631. 
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